Quantum Simulation of Molecular Collisions with Superconducting Qubits 
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We introduce a protocol for the fast simulation of n-dimensional quantum systems on n-qubit 
quantum computers with tunable couplings. A mapping is given between the control parameters 
of the quantum computer and the matrix elements of H s (t), an arbitrary, real, time-dependent 
n x n dimensional Hamiltonian that is simulated in the n-dimensional 'single excitation' subspace of 
the quantum computer. A time-dependent energy/time rescaling minimizes the simulation time on 
hardware having a fixed coherence time. We demonstrate how three tunably coupled phase qubits 
simulate a three-channel molecular collision using this protocol, then study the simulation's fidelity 
as a function of total simulation time. 



A quantum computer can significantly reduce the re- 
sources necessary to simulate quantum mechanical sys- 
tems pQ. Typically, quantum simulation algorithms al- 
gorithms construct the simulated system's time evolu- 
tion operator, energies and/or eigenstates from a uni- 
versal set of gates [2H5]. Alternatively, ultracold atoms, 
trapped ions, and liquid- state NMR have directly emu- 
lated the time evolution of certain other quantum sys- 
tems [9HTT]. Recent experimental progress suggests that 
quantum simulation will be one of the first practical ap- 
plications of quantum computation [7HT2]. 

In principle, an n-qubit quantum computer can store 
the state of any N = 2 n dimensional quantum system, an 
exponential reduction in the resources necessary to store 
quantum information on a classical computer. However, 
simulation may require ~ N 2 — 2 2n elementary gates per 
time step unless the simulated Hamiltonian has special 
properties, e.g. locality [2j [13]. Even for these special 
Hamiltonians, fully digital quantum simulation often re- 
quires an excessive number of gates for current quantum 
computing technology [4j |6] . 

In this Letter, we show that a subspace of a tunable 
n-qubit quantum computer can emulate an arbitrary n- 
dimensional quantum system, trading an exponential re- 
duction in resources for simulations of a wider variety 
of Hamiltonians. This subspace simulates other quan- 
tum systems very different from the computer itself in an 
amount of time that is independent of n. By comparison, 
classical simulation of an n-dimensional quantum system 
requires ~ n 3 elementary operations per time step. While 
the most efficient quantum simulation algorithms offer 
an exponential reduction in both qubits and elementary 
operations, they typically apply to specific, fundamental 
time-independent Hamiltonians, or those already simi- 
lar to that of the computer itself. We show that with a 
more modest polynomial reduction in resources, a sub- 
space of a tunable quantum computer can simulate any 
real, time-dependent Hamiltonian. 



We begin by outlining the theory behind our approach 
to simulation. First, we identify an n-dimensional invari- 
ant subspace suitable for quantum simulation. Then we 
define a time dependent energy /time rescaling that max- 
imizes the speed of the simulation within the constraints 
of the quantum computer. Finally, the control parame- 
ters of the quantum computer are given explicitly as a 
function of the matrix elements of H s (t). 

Our approach is tested by performing a simulation of 
a molecular collision with a circuit of tunably coupled 
Josephson phase qubits. Molecular collisions and elec- 
tronic structure calculations are widely studied as im- 
portant applications of quantum simulation techniques 
[BHS] • We show in detail how a superconducting circuit of 
three tunably coupled Josephson phase qubits simulates 
a three channel Na-He collision. Finally, we discuss the 
relationship between simulation fidelity and total simu- 
lation time for this particular example. 

An n- Dimensional Subspace of the full quantum com- 
puter's Hilbert space, can emulate another quantum 
system at all times only if it is invariant to the time evo- 
lution generated by the computer's Hamiltonian H qc (so 
that the subspace is well-isolated from the rest of T~L and 
evolves unitarily). We model i^ qc as 

ffqcto = E-^f + \y,9^)j,x ® ^ (i) 

i=l i^j 

where €i(t) are the uncoupled qubit energies, gij(t) = 
gji(t) are the pairwise qubit interaction strengths, J^ v 
gives the relative size of the <rf (g) a 1 - interaction, and 
/i, v G {0, x, n, z} are summed over. While €i(t) and gij(t) 
may in general be time-dependent, the time-independent 
structure of qubit interaction is specified by J^, a di- 
mensionless tensor that is typically fixed by a given ar- 
chitecture and is identical between each pair of qubits. 
In the weak coupling limit, \gij \\\J^ v \\/ei <C 1, subspaces 
of T~L are invariant to time evolution generated by i7 qc if 
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spanned by computational basis states having the same 
number of excited (tunable) qubits. The 'single excita- 
tion subspace', denoted as H n , is an n-dimensional in- 
variant subspace spanned by |z) = |00..01^..0 n ) for all 
i = 1, 2..., n. 

The control parameters t{ (t) and gij (t) directly control 
the Hamiltonian that % n simulates. We define H n as H qc 
projected into the single excitation subspace, 



H n (t) = PH qc (t)p1 



(2) 



where P is an nx2 n dimensional operator that projects T~L 
onto T-L n . Up to an additive energy shift, H n has matrix 
elements 



K j (t) 



- aY,k&9ik(t), i = j 
\9ij(t), i^j 



(3) 



with a = 2(J 20 + J zz ). We assume J xx + J yy ^ and 
normalize J^ v so that Jxx + Jyy = 1. In the weak coupling 
limit, H n is approximately invariant and generated by 



U n (t) EE PU qc (t)P^ 



J e -JifoH n {t>)dt> 



(4) 



where T is the time-ordering operator. generates U n 
exactly when no matrix elements of H qc mix l~i n with the 
rest of H (i.e. Jq x = Jo^ = J zx = = 0). The (n 2 + 
n)/2 parameters €i(t) and ^j(t) independently control 
each of the (n 2 + n)/2 matrix elements of the real H n 
and can therefore be used to simulate any arbitrary, real 
Hamiltonian in % n . 

While we can simulate H s in % n by choosing (t) and 
g%j(t) so that H n (t) = H s (t) for all £, a direct mapping 
between Hamiltonians limits the computer to simulating 
other quantum systems with similar energy scales over 
lengths of time within the computer's coherence time. 
Fortunately, simulation of H s only requires equality up 
to an overall phase between U n and the time evolution 
operator generated by H s : 



U(t) = Te~^^ Hs{t ' )dt ' 



(5) 



The time elapsed on the quantum computer, t qc (t), is a 
strictly increasing function of simulated time t, admit- 
ting a much less restrictive relationship between Hamil- 
tonians: 



H s (t) + c(t) = X(t)H n (t qc (t)). 



(6) 



c(t) is a time-dependent, additive energy shift giving the 
overall phase difference <j>(t) = ^ J t c(t')db' , and we have 
introduced a positive, time- dependent energy/time scal- 
ing 



X(t) dt qc /dt. 



(7) 



The energy /time scaling X(t) determines the speed of 
the simulation. By carefully minimizing A(t), we reduce 
the total simulation time and, consequently, the error 
due to decoherence. X(t) is bounded from below by ex- 
perimental constraints on the allowed values of control 
parameters €i(t) and gij(t) as well as their maximum 
rates of change. Suppose qubit interaction strengths can 
vary in a range g^t) G [-^ax,w], and the uncoupled 
qubit energies can vary in a range €i(t) G [e m i n , e max ]. For 
convenience, we define a simulated energy Ei(t) anal- 
ogous to €i(t) when diagonal contributions from qubit 
interactions are anticipated: 

Ei (t) ee Hf(t) + aY,H l J(t). (8) 

Using this definition together with equations (3) and (6), 
we relate the control parameters of the quantum com- 
puter to the simulated energies in H s (t): 



g l3 {t) = H^(t)/X(t) 
6i(t) = [E t (t)-c(t)]/X(t). 



(9) 



By choosing c(t) = E max (t) - A(t)e max where E max (t) is 
the largest value obtained by the Ej(t) at a particular t, 
we force each t{ to be as large as possible and therefore 
minimize leakage out of H n . 

Each of the computer's control parameters remains 
within its allowed range when \(t) is larger than (n 2 + 
n)/2 energy ratios at all times: 



\{t) > 



AEi(t)/A 

c max 



(10) 



where AEi(t) E max (t) - E^t) and Ae e max - e min . 
X(t) is also bounded by constraints on the speeds with 
which control parameters can change. Suppose v^(t qc ) = 
dei(t qc )/dt qc and v^(t qc ) = dgij(t qc ) /dt qc can never be 
larger in magnitude than v^ax an d v max respectively. 
Then for all £, 



v y > — 

max — ^2 



dH?(t) H?(t)d\ 



dt 



A dt 



(11) 



(and similarly for v^ ax ). 

To simulate H s (t) in < H n , we first choose X(t) as small 
as both inequalities (10) and (11) allow, guaranteeing a 
fast simulation within the experimental constraints of the 
quantum computer. We integrate over \(t) to calculate 
t qc as a function of t: 



f X(t')dt' + t qc {t { ). (12) 

With both X(t) and t qc (t) known, we can explicitly map 
the matrix elements of H s to the control parameters of 
the quantum computer: 



ei(t qc (*)) = 



AE z (t)/X(t) 



9ij (t qc (t)) = Hi*(t)/\(t). 



(13) 
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To demonstrate our theory in detail, we describe three 
Josephson phase qubits simulating a three-channel col- 
lision between a sodium and a helium atom. For three 
phase qubits with tunable inductive coupling, 

H qc (t) = + ^£ff«(*)*i®*i (14) 



= 1 



where 4>i is defined in terms of the matrix elements 
(fjk = (j\0i\k) of the local Josephson phase operator 



in the computational basis of the ith qubit: 



<P00 (fii ^ 
2^oi 



2(p 



Both the €i and the (fjk depend on the externally 
applied flux through the superconducting circuit. Ex- 
ternal flux bias is quantified by a dimensionless param- 
eter Si(t) = & x /&* where <£* is the qubit's critical flux 
bias, or alternatively, by the dimensionless well depth 
AU/ftwp [14]. We consider external bias values for which 
s e [.89, .90] and AU/hu p e [13.7,15.5]. In this range, 
Ae/h = 190MHz while ~ a\ + llof varies little. A 
tunable mutual inductance independently controls the 
couplings Qij(t) between each pair of qubits. We have 
assumed Josephson junction parameters Iq = 2.93 /iA, 
C = 1.52 pF, and L = 808 pH. 

An n-dimensional subspace can simulate a molecular 
collision only after we project the full, many-body Hamil- 
tonian of the interacting electrons and nuclei into an 
n-dimensional basis. We construct the collision Hamil- 
tonian from Born-Oppenheimer energies and nonadia- 
batic couplings calculated previously for three molecu- 
lar channels: Na(3s) + He(ls 2 ) [1 2 S+] and Na(3p) + 
He(ls 2 )[l 2 n ; 2 2 S+] [15 , labeled as |l) g , |2) g and |3) s 
respectively. The energies are stored for fixed values 
of the internuclear distance R, which we assume takes 
straight-line trajectories in a standard semiclassical ap- 
proximation: R(t) = \Jb 2 + v 2 t 2 where v is the incoming 
particle's velocity and b is the impact parameter of the 
collision. 

Figure 1 outlines our simulation protocol for H s (t) de- 
scribing a three-channel Na-He collision. The matrix el- 
ements of H s (t) are displayed in Fig. 1(a) for a given 
semiclassical trajectory R(t). Directly below, we plot 
the energy /time scaling parameter \(t) as a black curve 
enveloping the six energy ratios given in Eq. (10). A 
small X(t) speeds the quantum computer through times 
when the internuclear distance R is large, but as R de- 
creases (t —> 0), a relatively small g max value constrains 
the growing couplings. X(t) increases over two orders of 
magnitude, creating a highly nonlinear relationship be- 
tween t qc and t, as shown in Fig. 1(c). This effectively 
stretches the portion of the collision when internuclear 
distance is small over the entire simulation, as can be 
seen in the plot of the quantum computer's control pa- 
rameters as a function of t qc in Fig. 1(d). 
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FIG. 1: (color online) H s (t) describes a three chan- 
nel Na-He collision with b = 0.5 and v = 1.0. (a) 
Matrix elements of H s as a function of time in atomic 
units (Eh = 27.21 eV and the atomic unit of time is 
2.419 x 10 -8 ns). (b) The dimensionless time scaling 
parameter \(t) envelopes six energy ratios (AE3 = for 
all t). We assume g max /h = 2.0 MHz and Ae max //i = 
190 MHz. (c) Plot of t qc (t) for the case of t qc (t { ) = 0, 
t[ = —40 a.u.. (d) Control parameters that simulate 
H s (t) plotted as a function of t qc (63 = e max for all t qc ). 

To study the fidelity of the simulation, we compare the 
exact and simulated time evolution operators, U(t) and 
U n (t qc (t)) respectively, by plotting (in Fig. (2)) transi- 
tion probabilities out of |l): 



p u ( t ) = i(w)|i>i 2 



(16) 



Because the exact transition probabilities evolve differ- 
ently with t than the simulated evolve with t qc , we de- 
fine a time-dependent transition fidelity which accounts 
for time scaling, 

F(t) = \ s (l\uHt)U n (t qc (t))\l) n \ 2 , (17) 
and a time-dependent leakage out of H n , 

£(*) = £J±<i|M*qc(*))|l>J 2 (18) 
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FIG. 2: (a) (color online) Exact transition probabilities 
generated by H s (t) shown in Fig. 1(a). (b) Transition 
probabilities simulated with parameter profiles given in 
Fig. 1(d). Final simulation fidelity is 0.998. 

where J2 ± is the sum over all computational basis states 
orthogonal to 7~L n . In the upper part of Fig. (3), 
fidelity and leakage are plotted together for four different 
#max values. 



Minimizing g max \\Jn 



either by decreasing g n 



or by increasing e m i n , reduces leakage and thus improves 
simulation fidelity. In this example, we find simulation 
fidelity more sensitive to the cutoff in g max because leak- 
age is most prominent when the interatomic distances are 
small (t — >> 0) and the diabatic couplings between chan- 
nels are the dominant terms. By reducing g max we also 
increase X(t) and thus the total simulation time, as stud- 
ied in the lower plot of Fig. (3). To increase fidelity from 
.9990 to .9999 we need to increase the simulation time 
by a factor of ~ 3, a relationship that is independent 
of n. While not introducing specific models of decoher- 
ence, we note that high fidelity simulations are possible 
on superconducting qubits with coherence times around 
100 ns. 

When applied to molecular collisions, our approach to 
quantum simulation requires classical overhead to project 
the fundamental, time-independent, many-body Hamil- 
tonian into an independent, n-channel H s . The quanti- 
ties of physical interest, cross sections, are obtained by 
integrating the final transition probabilities over many 
semiclassical trajectories with different impact parame- 
ters, which requires no further classical overhead. A clas- 
sical simulation of transition probabilities requires ~ n 3 
elementary operations per time step for a single impact 
parameter, thus cross section calculations are computa- 
tionally intensive for large n. Alternatively, simulation 
time is independent of n using our protocol, so once the 
independent H s has been calculated, cross sections can 
be obtained quickly. 

In summary, we have presented a straightforward pro- 
tocol for quantum simulation that can be implemented 
with currently available superconducting quantum com- 
puting technology. While a promising application of 
quantum computation, current quantum simulation pro- 
tocols require a threshold number of gates and qubits 
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FIG. 3: (color online) (a) Fidelity and leakage as a 
function of simulation time for four different g max val- 
ues, with all other parameters the same as in Fig. 1. 
(b) Final simulation fidelity versus total simulation 
time for varying g max . The g max value is referenced by 
the shade of the data point. 

that prohibits fully digital quantum simulations from be- 
ing demonstrated on available quantum computers. How- 
ever, we have shown how quantum computers of only a 
few qubits can simulate arbitrary quantum systems ac- 
curately and quickly even before they reach the regime 
of fault tolerant quantum computation. 
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